library(tidyverse)
library(plotly)
library(DT)
library(limma)
library(clusterProfiler)
library(org.Hs.eg.db)
library(enrichplot)
library(WGCNA)
library(sva)
library(VennDiagram)
library(ggvenn)
LIMA
Rejection
rejection = readRDS("/Users/PengOlivia/Desktop/PROMAD-Download-Rejection-2024-08-17.Rds")
rejection = rejection$GSE98320
column_names <- colnames(rejection)
outcome = column_names[-c(1, length(column_names))]
rejection_vector <- ifelse(grepl("AR", outcome), 1, 0)
rejection_exprs_mat = as.matrix(rejection[, -c(1, 857)])
design <- model.matrix(~ rejection_vector)
limma.fit <- lmFit(rejection_exprs_mat, design)
limma.fit <- eBayes(limma.fit)
rejection_DEresults <- topTable(limma.fit, number = Inf, sort.by = "p")
rejection_DEresults$Gene <- rejection$X
DErejection <- rejection_DEresults[,c("Gene", setdiff(names(rejection_DEresults),"Gene"))]
datatable(DErejection)
Fibrosis
fibrosis = readRDS("/Users/PengOlivia/Desktop/PROMAD-Download-Fibrosis-2024-08-17 new.Rds")
fibrosis = fibrosis$GSE76882
column_names <- colnames(fibrosis)
outcome = column_names[-c(1, length(column_names))]
fibrosis_vector <- ifelse(grepl("Fibrosis", outcome), 1, 0)
fibrosis_exprs_mat = as.matrix(fibrosis[, -c(1, 143)])
design <- model.matrix(~ fibrosis_vector)
limma.fit <- lmFit(fibrosis_exprs_mat, design)
limma.fit <- eBayes(limma.fit)
fibrosis_DEresults <- topTable(limma.fit, number = Inf, sort.by = "p")
fibrosis_DEresults$Gene <- fibrosis$X
DEfibrosis <- fibrosis_DEresults[,c("Gene", setdiff(names(fibrosis_DEresults),"Gene"))]
datatable(DEfibrosis)
GSCE
Enrichment plot for
all rejection, fibrosis genes
#rejection
rejection_vector = DErejection$Gene
fibrosis_vector = DEfibrosis$Gene
biomarkers = intersect(rejection_vector, fibrosis_vector)
gsea_input = DErejection %>%
filter(Gene %in% biomarkers)
gene_list <- -log10(gsea_input$P.Value)
names(gene_list) <- biomarkers
gse <- gseGO(geneList=gene_list,
ont ="BP",
keyType = "SYMBOL",
OrgDb = org.Hs.eg.db)
datatable(gse@result)
path.id = gse@result$ID
gseaplot(gse, path.id[1])

#fibrosis
gsea_input = DEfibrosis %>%
filter(Gene %in% biomarkers)
gene_list <- -log10(gsea_input$P.Value)
names(gene_list) <- biomarkers
gse <- gseGO(geneList=gene_list,
ont ="BP",
keyType = "SYMBOL",
OrgDb = org.Hs.eg.db)
datatable(gse@result)
path.id = gse@result$ID
gseaplot(gse, path.id[1])

Venn diagram for top
5 overlapping genes
#venn diagram
rejection_vector_top_500 = DErejection[1:500, ]$Gene
fibrosis_vector_top_500 = DEfibrosis[1:500, ]$Gene
sets <- list(Set1 = rejection_vector_top_500, Set2 = fibrosis_vector_top_500)
ggvenn(
sets,
show_elements = F,
fill_color = c("red", "green", "blue")
)
